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ABSTRACT 

(~| ■ We use matter power spectra from cosmological hydrodynamic simulations to quantify the 

Q^l effect of baryon physics on the weak gravitational lensing shear signal. The simulations con- 

sider a number of processes, such as radiative cooling, star formation, supernovae and feed- 
back from active galactic nuclei (AGN). Van Daalen et al. (2011) used the same simulations 
^ , to show that baryon physics, in particular the strong feedback that is required to solve the 

. overcooling problem, modifies the matter power spectrum on scales relevant for cosmological 

weak lensing studies. As a result, the use of power spectra from dark matter simulations can 
lead to significant biases in the inferred cosmological parameters. We show that the typical 
biases are much larger than the precision with which future missions aim to constrain the dark 
energy equation of state, wq. For instance, the simulation with AGN feedback, which repro- 
duces X-ray and optical properties of groups of galaxies, gives rise to a ~ 40% bias in iuq. We 
■ also explore the effect of baryon physics on constraints on Jim, erg, the running of the spectral 

index, the mass of the neutrinos and models of warm dark matter We demonstrate that the 
modification of the power spectrum is dominated by groups and clusters of galaxies, the ef- 
fect of which can be modelled. We consider an approach based on the popular halo model and 
show that simple modifications can capture the main features of baryonic feedback. Despite 
its simplicity, we find that our model, when calibrated on the simulations, is able to reduce 
the bias in wq to a level comparable to the size of the statistical uncertainties for a Euclid-like 
mission. While observations of the gas and stellar fractions as a function of halo mass can be 
kJ ' used to calibrate the model, hydrodynamic simulations will likely still be needed to extend 

'^-^ ' the observed scaling relations down to halo masses of 10^^ h^^ Mq. 
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1 INTRODUCTION 

The discovery that the expansion of the Universe is accelerat- 
ing is arguably one of the most significant discoveries of modem 
cosmology. We lack, however, a theoretical framework to explain 
the observations. For instance, current measurements can be ex- 
plained by a cosmological constant or dynamic mechanisms (e.g., 
quintessence). The field driving the acceleration has been dubbed 
"dark energy", which can be parametrised by an equation-of- state 
liiQ (and additional evolving parameters). Modifications of the laws 
of gravity on cosmological scales have been considered as well (for 
a recent review, see for example Uzan 2007). 
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Observations of type la supernovae provided the first 
evidence for the accelerated expansion ( IPerlmutter et al. 19991 
IRiess et al. 20071 lAmanuUah et al. 2010b . Further support has 
come from observations of the cosmic microwave background 
(e.g., Komatsu et al. 2011), baryon acoustic oscillations (e.g., 
Eisenstein et al. 2005; Percival et al. 2010) and gas fractions of 
massive clusters (e.g., Allen et al. 2008). These methods allow us 
to constrain the properties of dark energy solely by studying the 
expansion history of the Universe and do not distinguish between 
modified theories of gravity that do not alter the expansion history. 

Studies of the growth of structure are sensitive to both dark 
energy and modifications of gravity, and are therefore of great in- 
terest. The interpretation of the observations is complicated by the 
fact that most of the matter is in the form of dark matter. Fortu- 
nately, inhomogeneities of the matter distribution cause the light of 
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background sources to be differentially deflected. This leads to co- 
herent alignments in the shapes of these galaxies, which can be re- 
lated directly to the power spectrum of matter density fluctuations. 
This phenomenon is referred to as "weak gravitational lensing" and 
when applied to the study of large-scale structure, it is also known 
as "cosmic shear". 

The cosmic shear signal has been measured using large 
ground-based surveys (e.g., Hoekstra et al. 2006; Fu et al. 2008) 
and from space (e.g., Massey et al. 2007; Schrabback et al. 2010). 
As discussed by the Dark Energy Task Force (DETF; Albrecht et al. 
2006), it is one of the most promising probes of dark energy. Space- 
based missions are currently planned with the aim to constrain the 
dark energy equation of state with a relative precision of ~ 1%. The 
interpretation of the signal will, however, require exquisite knowl- 
edge of the power spectrum, well into the non-linear regime. 

In order to exploit the statistical power of these future data 
sets, we need to be able to model the power spectrum with percent 
accuracy. To date, the interpretation of cosmic shear measurements 
rely on results based on N-body simulations that contain only colli- 
sionless (dark matter) particles. Under the assumption that the evo- 
lution of the large-scale structure is mainly described by the evo- 
lution of the dark matter distribution, this should lead to fairly ac- 
curate results. Recently, Heitmann at al. (2010) developed an 'em- 
ulator' which is able to provide dark matter only models with 1% 
accuracy on scales larger than a megaparsec. 

However, a complete description of the matter fluctuations 
should also include baryon physics. Bernstein (2009) suggests to 
account for the effect of baryons on two-point weak lensing mea- 
surements by adding an unknown component to the power spec- 
trum which is expressed using Legendre polynomials. The nuisance 
of the baryons is then included in the figure of merit (Albrecht et al. 
2006) by marginalising over this component. Similarly, Kitching & 
Taylor (2010) propose a path-marginalisation technique and claim 
that this procedure reduces the figure of merit by 10%. In con- 
trast, Zentner, Rudd & Hu (2008) suggest that constraining bary- 
onic feedback and cosmology at the same time would degrade the 
figure of merit by about 30% — 40%. From those studies it emerges 
that self-calibrating techniques are very promising although the 
strength of the result relies on prior assumptions; because of that 
it seems clear that the better we understand how various mecha- 
nisms affect the density distribution of the baryons, the better we 
can choose priors in the self-calibration procedures. For this reason 
it is important to study the effect of baryonic feedback directly us- 
ing hydrodynamic simulations and this is the approach we adopt in 
this paper. 

Hydrodynamic simulations are much more expensive than col- 
lisionless N-body simulations, and it has only recently become 
possible to simulate cosmological volumes with reasonable res- 
olution. These simulations aim mainly to reproduce the proper- 
ties of galaxies and the intergalactic medium, and to do that they 
need to account for processes such as radiative cooling, star for- 
mation, stellar mass loss, chemical enrichment and outflows driven 
by supemovae and AGN. To complicate the situation further, it is 
not yet clear which processes one needs to include in the simu- 
lations to mimic the observed Universe and the simulations still 
lack the resolution to model all relevant processes from first prin- 
ciples. The Overwhelmingly Large Simulations (OWLS) project 
( [Schaye et al. 2010l l represents an effort in that direction. The main 
goal of the OWLS project is to study how the properties of large- 
scale structures, galaxies and the intergalactic medium change in 
different scenarios, by producing a large number of simulations 
and by varying parameters independently. Comparisons with ob- 



servations allow one to establish which mechanisms are required to 
better describe our Universe. One of the main results of the OWLS 
project is that the simple insertion of baryons, star formation and 
feedback from supemovae, cannot reproduce the characteristics of 
groups of galaxies. Accounting for AGN feedback, however, re- 
sults in good agreement with both X-ray and optical observations 
(McCarthy et al. 2010). 

Importantly for cosmic shear studies, van Daalen et al. (201 1) 
used the OWLS runs to show that the resulting power spectrum of 
matter fluctuations depends strongly on the adopted baryon physics 
out to surprisingly large scales. In particular, their model for AGN 
feedback leads to significant changes compared to scenarios with 
only star formation and relatively inefficient feedback from super- 
novae. In this latter scenario, cooling baryons change the profiles 
of the haloes, generating more compact cores and increasing power 
on small scales (e.g., Jing et al. 2006; Rudd et al. 2008, Guillet et al. 
2010; Casarini et al. 201 1). To prevent the formation of overly lumi- 
nous galaxies and to reproduce the X-ray properties of groups, ex- 
tremely efficient feedback is required. AGN feedback can prevent 
the overcooling of the baryonic component by ejecting large quan- 
tities of gas at high redshift, when the supermassive black holes 
were growing rapidly (McCarthy et al. 2010; 201 1). Consequently, 
AGN feedback can modify the power spectrum in a dramatic way, 
decreasing the scale below which cooling increases the power, and 
significantly reducing the power up to very large scales. Levine & 
Gnedin (2006) used a simple toy model to manually include AGN 
outflows in collisionless simulations and to show that AGN can po- 
tentially drastically change the power spectrum, although they were 
unable to predict the sign and the magnitude of the effect. More 
recently, van Daalen et al. (2011) predicted, using a fully hydro- 
dynamic simulation that reproduces optical and X-ray observations 
of groups of galaxies, a ~ 10% decrease in total power relative to 
a dark matter only simulation on scales of 0.1 — 10 h^^ Mpc at 

2 = 0. 

Ignoring the effects of baryon physics may therefore lead to 
biases in the constraints of cosmological parameters inferred from 
cosmic shear studies. In this paper we use the power spectra tabu- 
lated by van Daalen et al. (201 1) to evaluate the effects of a number 
of feedback scenarios. In addition, we provide a simple, yet rather 
effective approach to reduce the bias of cosmological parameters. 
We show that this kind of approach will be crucial for upcoming 
surveys because the bias will be much larger than the statistical un- 
certainties, if left uncorrected. 

This paper is organised as follows. In section |2] we describe 
the set of simulations we use for this work. In section|3]we examine 
how cosmic shear statistics are affected. In section |4] we quantify 
the biases introduced when ignoring the effects of baryon physics. 
In section [5] we show that it is possible to recover realistic power 
spectra and to remove (most of) the biases by making use of simple 
models that are calibrated to match the star and gas fractions of the 
simulated groups. Cosmic shear cannot only constrain the proper- 
ties of the dark energy, but can also be used to study the running 
of the spectral index, warm dark matter and massive neutrinos. In 
section [6] we evaluate the consequences of baryonic feedback on 
the ability of cosmic shear to constrain those effects. Finally, we 
conclude in sectionl?] 



2 SIMULATIONS 

To examine the physics that drives the formation of galaxies and the 
evolution of the intergalactic medium, the Overwhelmingly Large 
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Simulations (OWLS) project includes over 50 large, cosmological, 
hydrodynamic simulations run with a modified version of the SPH 
code Gadget (last described in Springel 2005). A range of physical 
processes was considered, as well as a range of model parameters. 
In this paper we use a subset of OWLS for which we will give a 
brief description, inviting the reader to find more details in the pa- 
pers where the simulations are presented. Note that all simulations 
have been performed with the same initial condition^ The sim- 
ulations considered here are (following the naming convention of 
Schaye et al. 2010): 

• DMONLY: a dark matter only simulation, of the kind com- 
monly used to compute the non-linear power spectrum which is 
needed in weak lensing studies. It is therefore the reference to 
which we compare the other simulations. 

• REF: although it is not the reference simulation for the study 
presented here, this simulation includes most of the mechanisms 
which are thought to be important for the star formation history (see 
Schaye et al. 2010 for a detailed discussion) but not AGN feedback. 
The implementation of radiative cooling, star formation, super- 
novae driven winds, and stellar evolution and mass loss have been 
described in Wiersma, Schaye & Smith (2009), Schaye & Dalla 
Vecchia (2008), Dalla Vecchia & Schaye (2008), and Wiersma et 
al. (2009), respectively. This simulation represents a standard sce- 
nario assumed in cosmological hydrodynamic simulations. 

• DBLIMFV1618: this simulation has been produced using the 
same mechanisms as REF. The only difference between the two 
simulations is that in this simulation the stellar initial mass func- 
tion (IMF) was modified to produce more massive stars when the 
pressure of the gas is high, i.e. in starburst galaxies and close to 
galactic centres. This is obtained by switching from the Chabrier 
(2003) IMF assumed in the REF model to a Baugh et al. (2005) 
IMF in those regions. There are both observational and and the- 
oretical arguments to support a top-heavy IMF in those extreme 
conditions (e.g. Padoan et al. 1997; Baugh et al. 2005; Klessen et 
al. 2007; Maness et al. 2007; Dabringhausen et al. 2009; Bartko et 
al. 2010; Weidner et al. 2010). The IMF change causes the number 
of supemovae and the effect of stellar winds to increase resulting in 
a suppression of the SFR at smaller redshifts. However, this mech- 
anism alone is not able to reproduce the observed SFR (see Schaye 
et al. 2010). 

• AGN: a hydrodynamic simulation which differs from REF 
only by the inclusion of AGN. The AGN feedback has been mod- 
elled following Booth & Schaye (2009). In this approach AGN 
transfer energy to the neighbouring gas, heating it up and driv- 
ing supersonic outflows which are able to displace a large quan- 
tity of baryons far from the AGN itself. Among the three simula- 
tions considered here, it is arguably the most realistic, as it is able 
to reproduce the gas density, temperature, entropy, and metallic- 
ity profiles inferred from X-ray observations, as well as the stel- 
lar masses, star formation rates, and stellar age distributions in- 
ferred from optical observations of low-redshift groups of galaxies 
( [McCarthy et al. 20I0| l. 

To forecast the cosmic shear signal for the four different sce- 
narios, we make use of the results of van Daalen et al. (201 1), who 
tabulated the power spectra of matter fluctuations P{k, z) in red- 
shift slices over the redshift range ^ z ^ 6 for a number of 



The cosmology used to realise the simulations is the best-fit to 
the WMAP3 data I Spergel et al. 2007i: {^m, f^b, f^A, o-8, "s, 'i} = 
{0.238, 0.0418, 0.762, 0.74, 0.951, 0.73}. 
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Figure 1. Ratio between the power spectrum of matter fluctuations mea- 
sured from the simulations with baryons and the one measured from the 
DMONLY simulation. The ratio for the REF simulation is shown in green, 
the one for the AGN simulation is shown in blue, and the one for the 
DBLIMFV1618 model is shown in pink. Since the simulations have been 
carried out using the same initial conditions, deviations of the ratio from 
unity are due to the differences in baryon physics. 



OWLS runs. A detailed discussion of the procedure to compute the 
power spectrum and its accuracy can be found in van Daalen et al. 
(2011). Their convergence tests and noise estimations suggest that 
the power spectra estimated from OWLS are reliable up to at least 
fc ~ 10 /iMpc~^ over the range of redshifts we are interested in 
(i.e. z < 1, as the lensing signal is most sensitive to structures that 
are halfway between the observer and the source). At small fc, the 
estimate of the power spectrum is affected by the finite size of the 
simulation box (100 Mpc on a side). This is not a concern, be- 
cause on these scales baryonic effects are very small and density 
fluctuations are in the linear regime, so that we can compute the 
power spectrum from theory instead. 

Figure[T]shows the power spectrum measured for each simula- 
tion in three redshift bins normalised by the power spectrum of the 
dark matter simulation (DMONLY) at the same redshift. In the REF 
scenario (green), the presence of the baryons slightly suppresses the 
power spectrum at intermediate scales, due to the pressure of the 
gas. At smaller scales where baryons cool, the power spectrum is 
enhanced as the baryons fall into the potential wells. For this model, 
only the small scales are affected in an almost redshift indepen- 
dent way. The effect of baryon physics is more pronounced for the 
DBLIMFV1618 model, and depends on redshift. The AGN model 
leads to the largest difference compared to the DMONLY simula- 
tion. The amplitude of the power spectrum is strongly reduced on 
scales of ~ 1 — 10 Mpc and the effect increases as the redshift 
decreases; this is in agreement with the results by McCarthy et al. 
(2011) who showed that because AGN remove low-entropy gas at 
early stages (2 < 2 < 4), the high-entropy gas left in the haloes 
does not cool down and form stars and the suppression of power 
becomes more and more accentuated at small scales. 

The latter two scenarios are qualitatively similar, although 
the mechanisms are different: in the DBLIMFV1618 simulation 
baryons are removed due to the enhanced supernova feedback, 
whereas in the AGN scenario they are removed mostly by AGN 
feedback, at least for the more massive and thus strongly clustered 
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Figure 2. Top panel: relation between the wave number fc = s / f (wmux) 
and angular wave number s for vaiious source redshifts. This relation 
shows, for a given source redshift. which k contributes most to the con- 
vergence power spectrum at a given s. For both s and k we show the 
corresponding real-space conjugate variables. Bottom panel: typical Pk{s) 
signal-to-noise ratio as a function of the wave number s. We show results 
for the same set of source redshifts as were shown in the upper panel. 

haloes. Thus, the fraction of baryons which is removed is different, 
as is the rate at which they are removed. At very small scales cool- 
ing still enhances structure formation, but the physical scale below 
which this occurs is smaller than in the REF simulation. 

Although we will focus on a set of simulations which have 
been produced using the WMAP3 cosmology, van Daalen et al. 
(2011) compared a dark matter only simulation and an AGN sim- 
ulation realised with the best-fit WMAP7 cosmology (Komatsu et 
al. 2011) and found that the difference is the same as between the 
AGN and DMONLY simulations used here. This implies that the 
relative effect of baryonic feedback does not change significantly 
with the cosmology, and thus that our conclusions are not restricted 
to a specific set of cosmological parameters. 



3 EFFECT OF BARYONS ON TWO-POINT SHEAR 
STATISTICS 

In this section we briefly introduce the basics of weak gravitational 
lensing by large-scale structure and how baryon physics affects the 
interpretation of the measurements. For a more extensive review of 
cosmic shear, see for example Hoekstra & Jain (2008) and Munshi 
et al. (2008). A thorough discussion of the theory of weak lensing 
and its applications is given in Bartelmann & Schneider (2001). 

Massive structures along the line of sight deflect photons 
emitted by distant galaxies. Provided the source is small, the ef- 
fect is a remapping of the source's surface brightness distribution: 
the source is both (de)magnified and sheared. In the weak lens- 
ing regime, the convergence n gives the magnification (increase 
in size) of an image and the shear 7 gives the ellipticity induced 
on an initially circular image. Under the assumption that galaxies 



are randomly oriented in the absence of lensing, the strength of 
the tidal gravitational field can be inferred from the measured el- 
lipticities of an ensemble of sources. The resulting complex shear 
7 = 71 + 172 is a spin-2 pseudo vector, which can also be written 
as 7 = 7 exp(2ici), where a is the position angle of the shear. 

If redshift information is available for the sources, one can 
compute the two-point shear coiTelation functions for galaxies in 
the i*'' and j"* redshift bins, which are defined as 

= (7.,t(0i)7.,t(02) +7.,x(0i)7.,x(02)> (1) 

C'W = (7.,t(0i)7j,t(02) - 7.,x(0i)7.,x(02)> (2) 

where 6 = |0i — 62! and we have defined the tangential and 
cross components of the shear as: 7t = — 7?.e(7 exp( — 2i(/))) and 
7x = — Xm(7 exp( — 2i0)), with tj) the polar angle of the sepa- 
ration vector 0. Note that the ensemble average depends only on 
the angular distance between the galaxies. The measurement of the 
redshift dependence of the cosmic shear signal greatly improves the 
constraints on cosmological parameters from weak lensing and is 
a key goal of current and future surveys. It is often refeiTed to as 
weak lensing tomography (e.g. Hu 1999; Hu 2002), because it al- 
lows us to study the matter distribution in "slices". The correlation 
functions can be measured from a catalogue of galaxy shapes and 
they are related to the convergence cross-power spectrum: 

1 r°° 

CVW = ^/, dsVJo/4(s'e)P:^(.'), (3) 

where Jo, J4 are the O"^ and the 4*'' order Bessel functions of 
the first 8kind. The convergence power spectrum P^-' (s) is re- 
lated to the power spectrum of matter fluctuations P{k, w) through 
dKaiser 1998llSchneider et al. 1998t : 

4c4 Jg a^(w) V/k(w) / 

with 

gi(w)^ dwp^(w) — — , (5) 

where /k (w) is the comoving angular distance, w is the radial co- 
moving coordinate, wu is the radial comoving coordinate of the 
horizon; Hq, Sim and a{w) are the Hubble constant, the matter 
density parameter and the scale factor, respectively. The projected 
power spectrum depends on Pi{w), the radial distribution of the 
sources in the i"* redshift bin. Since s is the Fourier-conjugate of 
the angle 6, we can relate an angular scale to it through s = 2n/6. 

Various two-point shear statistics that have been employed 
in the literature simply correspond to different filters of Pk{s) in 
Eqn.|4] Therefore, the detailed effect of baryonic feedback will de- 
pend somewhat on the statistic that is used, but this does not signif- 
icantly affect our conclusions which are derived using Pk (s). Note 
that in the case of a full sky survey one could determine the C{1) 
from a decomposition in spherical harmonics. For small angular 
separations s I and /(/ + 1)C(/) s^P„(s). 

3.1 Relevant scales 

Cosmic shear results are typically presented in terms of angular 
scales, rather than the physical scales shown in Fig. [T] We there- 
fore start by examining which angular scales are affected by baryon 
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Having established that cosmic shear studies are sensitive to 
the scales where baryon physics modifies the power spectrum, we 
now want to quantify how various scenarios change the two-point 
shear statistics. For that we adopt a source redshift distribution that 
is representative of the CFHTL S-Wide ([Benjamin et al. 2007) and 
a fair approximation for Euclid ( jLaurejis et al. 2009) . We adopt the 
following parametrisation: 



p{z) 



(6) 



with a = 0.836, /3 = 3.425, and 20 = 1.171. We divide the 
source galaxies in three tomographic bins with limits [0, 0.6, 1.2, 
3.4], which yields six cross-power spectra. 

The top panel of Figure|3]shows the value of ^+ [6) measured 
for the various feedback scenarios, normalised by the results for 
DMONLY. The effect of baryons is small and limited to very small 
scales for the REF scenario. However, for DBLIMFV1618, and in 
particular for the AGN model, the difference with the DMONLY 
result is large and increases when the redshift of the sources de- 
creases. The redshift dependence is the result of two effects. The 
first is a geometric one: when the redshift of the sources decreases, 
the physical scales probed by the lensing signal become smaller 
(see Figure|2j. The second reason is the suppression of the ampli- 
tude of the power spectrum due to feedback, which becomes larger 
at late times (see Figure [T). The bottom panel of Figure |3] shows 
the value of ^_ [6) measured for the various feedback scenarios, 
normalised by the results for DMONLY. Notice that the bias for 
5_ is more pronounced out to larger scales. This is because ^_ is 
much more sensitive to small-scale structures (i.e. to the shape of 
the power spectrum Pk{s) for large s). 



Figure 3. Top panel: ratio of the correlation function ^+{6) for 
REF/DMONLY (green), DBLIMFV1618/DMONLY (pink) and 
AGN/DMONLY (blue) . The notation hiriij indicates the correlation 
of sources from redshift bin i with sources from redshift bin j. Here, we 
show only results from the bins with i = j. Bottom panel: same as the 
upper panel but for the correlation function f _ (8). 



physics. To do so, we consider sources at a single redshift. In this 
case one can show that the value of Pk{s) depends mostly on the 
density fluctuations with comoving wave numbers ~ s//A'(ifmax) 
with fK{wms.^) maximising the ratio ^^''^f~^^^^^'' ■ In the top 
panel of Figure |2l we show, for various source redshifts Zs, the re- 
lation between the angular wave number s and the wave number 
s//K(wmax) uslug the adopted WMAP3 cosmology. It shows, for 
example, that measuring the power spectrum P„ (s) at s ~ 1 x 10* 
of galaxies with redshifts ~ 0.8, probes density fluctuations at 
scales k ~ lOh Mpc~^, where baryon physics is important. 

However, one might wonder if the signal at arcminute scales 
is statistically important. To examine this, the bottom panel of Fig- 
ure |2] shows a typical signal-to-noise ratio of Pk{s). The signal 
has been computed assuming a WMAP3 cosmology. The noise ac- 
counts for sampling and statistical noise, assuming a WMAP3 cos- 
mology and a survey area A = 20000 deg^, a number density 
of galaxies of n = 30 gal/arcmin'^ all placed at the same red- 
shift Zs and with intrinsic ellipticity dispersion ae = 0.33 (see sec- 
tion |4] for more details on the noise computation). As one can see, 
the signal-to-noise ratio peaks at scales between 2 and 10 arcmin, 
where baryon physics is important. 



3.2 Effect on cosmological parameter estimation 

It is clear from Figure [T] that the change in the power spectrum 
is large in the case of the AGN and DBLIMFV1618 scenarios. 
The modification is, however, scale-dependent, which may help to 
ameliorate the problem, since this cannot be reproduced by vary- 
ing cosmological parameters which predominantly affect the over- 
all amplitude of the weak lensing power spectrum. In other words, 
it might be possible to separate the effects of baryonic feedback, or 
at least to identify them: the inferred values for cosmological pa- 
rameters from weak lensing statistics are scale-dependent for the 
AGN and DBLIMFVI6I8 scenarios. 

We first investigate the effect on the recovered value of as, the 
rms fluctuation of matter in spheres of size 8 h^^ Mpc. A compli- 
cation to our analysis is the limited accuracy of the prescriptions 
for the non-linear power spectrum, be it Peacock & Dodds (1996) 
or the halofit approach (Smith et al. 2003) used here. We therefore 
cannot predict 5+,dmonly(6', Zs) directly, but the procedure out- 
lined below is accurate as the predictions should have the correct 
scaling as a function of as- For the various feedback models we 
first define the ratio 

Z? fa ^ — C+,hydro(^, Zs) 

4+,DMONLY(t7, Zs ) 

as a function of source redshift Zs and angular scale 9. Here 
C+,hydro(6', 2s) is the Correlation function measured for REF, 
DBLIMFVI6I8 or AGN, whereas C+,dmonly(6', 2s) is the 
DMONLY correlation function. We use the halofit prescription 
(Smith et al. 2003) to compute ^+,haiofit{d, Zs; as), keeping all 
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Figure 4. Top (bottom) panels show the deviation of the inferred trg (wq) from the true reference value a-8,rcf = 0.74 (uio,rof = —1) as a function of source 
redshift, when the amplitude of the ellipticity coiTelation function ^-|- (6) is used to estimate the cosmological parameter of interest (while the other parameters 
are kept at their reference values) and when we use halofit models (see text for details). The deviation depends on the angular scales that is used and is smaller 
for larger scales. The left panels show the results for the REF scenario, the middle panels for the DBLIMFV1618 and the right panels for the AGN scenario, 
which results in the largest biases. 



other cosmological parameters fixed to the reference values. We 
define the ratio 

P (a ^ .„\ C+,haloflt(6',2:s;0-8) 

-«+, halofit (b', 2s, O-gj = -J-Z r-, (S) 

? + , ref, halofit (f, Zs) 

where ^+,rcf,haiofit is the halofit prediction for the reference cos- 
mology. For each Zs and 9 we find the value of erg for which 
R+,ha.iofLt{6, Za) = -R+,hydro(f , Zs). Hence, we compute by how 
much the value of needs to change from the reference value if 
one ignores feedback processes and instead interprets the measure- 
ment of ^+ (S, ) in a dark matter only framework. 

The top panels of Figure |4] show the resulting value Acrg = 
as — (Tg,ref for which 7?haioflt = Rhydio 38 a functlon of Zs for vari- 
ous values of As we anticipated, the infeiTed values of ug depend 
on angular scale, but the effect of baryon physics is modest, partic- 
ularly for the REF simulation. For the AGN and DBLIMFV1618 
models the effect is qualitatively the same and it is still only a few 
percent. 

In models with a constant dark energy equation of state P = 
Wop, the change of wq mainly leads to a change in the amplitude 
of the weak lensing power spectrum. It is therefore interesting to 
repeat the same analysis for models with dark energy. Note that if 
we vary the value of wo, both the expansion history and the history 
of structure formation change. Because we normalise the amplitude 
of the fluctuations at the present time, decreasing the value of wo 
results in fluctuations that are larger at earlier times. On the other 
hand, a very negative equation of state increases the expansion rate. 



Thus, the overall amplitude of the two-point shear statistics at a 
given scale and redshift depends on which of the two competing 
effects is dominant. 

The bottom panels of Figure |4] show the values Awo = 
^"0 — ifo,rcf as function of the angular scale and source redshift. 
Compared to the bias in ag, the change in the inferred value of too 
is more dramatic, reaching 20% for the AGN scenario, 10% for 
the DBLIMFV1618 scenario, and even for the REF scenario wo 
is few percent too high. Note that we did not consider a redshift 
dependent equation of state w{z); also in that case we expect the 
estimated value of w to change as a function of the angular scale, 
leading to similar conclusions. These findings suggest that baryon 
feedback can lead to significant biases in cosmological parameters 
in the case of future cosmic shear studies. In the next section we 
will quantify this in more detail. 



4 EFFECT ON LIKELIHOOD RESULTS 

In this section we quantify the bias in the inferred values for wo and 
erg if one interprets weak lensing measurements using dark matter 
derived models. To do so, we perform a likelihood analysis, where 
we define the posterior probability distribution as: 

P(p|d) cx P(p)z:(d|p) (9) 

where p and d are vectors of parameters and observed data, respec- 
tively, and the likelihood £(d|p) is given by: 
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Figure 5. Left panel: f2m — likelihood contours for a CFHTLS-like survey. Middle panel: — wq likelihood contours for a CFHTLS-like survey. Right 
panel: posterior probability distribution for p(«>o) obtained for erg = 0.74 and Hni = 0.238. SoHd contours mark the 68% confidence regions. The shifts 
relative to the DMONLY case indicate the presence of bias due to baryonic effects. The biases are largest for the AGN scenario. 
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Figure 6. Left panel: joint constraints on erg — Qm for a Euclid-like survey. Solid, dotted and dashed lines mark the 68%, 95% and 99% confidence regions 
respectively. Middle panel: joint constraints on Hm — uiq. Right panel: posterior probability distribution for piwo) marginalised over erg and Qni-The shifts 
relative to the DMONLY case indicate the presence of bias due to baryonic effects. The biases are largest for the AGN scenario. Note that for the right panel 
the sign of the shift differs from what was found in Figure [s] because we now marginalise over the other parameters, rather than keeping them fixed. 



£(d|p) oc exp [-l/2(m(p) - d)C"^(m(p) - d)^] . (10) 

Here m(p) is a model and C is the covariance matrix. We chose 
the data vector d to be Pui^s) sampled at twenty scales between 
s = 10 and s = 6000 for the six cross spectra introduced in 
section [5] The total data vector thus contains 120 measurements. 
We compute cosmological models Pk{s) at the same scales using 
the halofit prescription for the non-linear power spectrum. All the 
models use the transfer function by Eisenstein & Hu (1998) which 
accounts for the effect that baryons have on the power spectrum 
during the radiation-dominated epoch. In this way the linear power 
spectrum of our models is similar (Eisenstein & Hu is just an ap- 
proximation) to the one used to establish the initial conditions of 
the OWLS simulations. 

Throughout the paper we use flat priors — [0.160, 0.316] 
and (Tg — [0.65,0.83], corresponding to the ±3(t error-bars for 



WMAP7. We also adopt a flat prior of wo = [-2.00, -0.6]. We 
emphasise that uninformative priors are needed to study the bias, 
because (unbiased) information from external data may force the 
recovered values towards their unbiased values. 

Following Takada & Jain (2004), we compute the covariance 
matrix in the Gaussian approximation. Note that the Gaussian ap- 
proximation breaks down at small scales. This is because the non- 
linear evolution of the density field causes the modes to mix, result- 
ing in non-zero off-diagonal terms. As a consequence, the value of 
the variance is increased (Semboloni et al. 2007; Takada & Jain 
2009; Pielorz et al. 2010). Since we do not want to overestimate 
the impact of baryon physics by underestimating the error bars at 
small scales, we only consider modes s ^ 6000, which correspond 
to angular scales larger than 9 ~ 3.5arcmin. 

As before, we do not use the data vector dhydro, which is mea- 
sured from the simulations. Instead we construct a new data vector 
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dfit: 



J dhaloflt J 

dflt = -2 X d 



hydro ; 



(11) 



100 



6min(arcmin) 



10 



where dhydro corresponds to Pk{s) from the REF, AGN, or 
DBLIMFV1618 simulation, and where dDMONLY is computed us- 
ing the DMONLY simulation. The halofit approach is used to cal- 
culate dhaloflt, which quantifies the dependence on the cosmology. 
If the halofit model described the power spectrum of the DMONLY 
simulation perfectly then the ratio in Eq.[TT]would be unity for any 
component of the data vector and d^t would be merely dhydro- 
This is, however, not the case because the halofit model has lim- 
ited accuracy; moreover the set of simulations we use is relatively 
small so the measured power spectrum is affected by sampling vari- 
ance. We are only interested in the comparison between dark matter 
only simulations and simulations with baryonic feedback. By using 
Equation [TT] we can compute the bias in cosmological parameters 
that one would obtain by neglecting the existence of baryons, while 
minimising the limitations of the halofit prescription and sampling 
variance due to the finite size of the simulations. 

We first consider a survey with area A = 200 deg^, galaxy 
number density n = 15 gal/arcmin^ and intrinsic ellipticity dis- 
persion of galaxies ae = 0.44, which roughly corresponds to the 
CFHTLS-Wide survey (e.g., Hoekstra et al. 2006; Fu et al. 2008). 
The joint constraints on Qm and as are presented in the left panel 
of Figure|5] whereas the constraints on f2,n and wq are displayed in 
the middle panel. Note that we keep the other parameter(s) fixed, 
rather than marginalising over them, because the statistical power 
of the CFHTLS is too limited to constrain them. The shifts in the 
likelihoods, compared to the DMONLY contours, indicate that the 
cosmological parameters are biased if a dark matter power spec- 
trum is used to interpret the data. In agreement with our earlier 
findings, the shifts are small for the REF model and larger for the 
other two scenarios. 

Based on the results of the previous section, we are most con- 
cerned about a bias in wq. The right panel of Figure |5] shows the 
posterior probability for wq, when the other parameters are kept 
fixed at their fiducial values. This result suggests that a dark matter 
only modelling of the CFHTLS signal will likely predict a value of 
Wo which is about 10% less negative than the real one. However, 
the DMONLY peak lies at the edge of the 68% confidence region 
of the AGN scenario, thus the bias is similar to the statistical un- 
certainty. Note that since the CFHTLS-Wide has limited statistical 
power due to its modest survey area, we fixed the values of as and 
Qui to constrain too, which affects the bias. We will see below that 
the bias is different (and increases) if we use a flat prior. 

To reach a precision of a few percent, much larger projects 
are being planned. Of these, we focus on the space-based Eu- 
clid mission ( [Refregier et al. 2010) , which aims to survey an area 
A = 20000 deg with a number density of galaxies of n = 
30 gal/arcmin^ and intrinsic ellipticity dispersion tJe — 0.33. This 
leads to much smaller statistical errors on the cosmological param- 
eters. As shown in Figure|6] in this case the effect of baryons leads 
to significantly biased results. The left panel shows the constraints 
on fini and as, while marginalising over the value of wo- The like- 
lihood contours for both the AGN and DBLIMFV1618 models are 
shifted towards lower values for fim, and an only slightly higher 
as- The consequences for the constraints on wq are even more dra- 
matic, as is clear from the middle panel of Fig. |6] The bias in wq 
reaches almost 40% for the AGN scenario, as is also evident from 
the marginalised probability distribution shown in the right panel. 
The large bias is the consequence of the large modification of the 
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Figure 7. Difference Awq between the best fit value of wg and the true ref- 
erence value «)o = — 1 as a function of the maximum angular wavenumber 
•5max (or minimum scale 6'njin) that is included in the likelihood analy- 
sis. The error bars represent the resulting la uncertainties on wq. Avoiding 
high wave numbers allows one to reduce the bias affecting the cosmological 
parameters at the cost of an increase of the statistical errors. 



power spectrum at intermediate scales. Note that the shift in the 
value of Wq in the right panel is in the direction opposite to what 
was found in Figure|4] because here we marginalised over the other 
parameters whereas they were kept fixed in Figures 2]and|5] 

Since the bias arises due to differences in the power spectrum 
at relatively small angular scales 6, or large s, it is interesting to 
examine whether the bias can be reduced by leaving the large wave 
numbers out of the likelihood analysis. Figure |7] shows the maxi- 
mum likelihood value of wq as a function of the maximum angular 
wave number Smax, while marginalising over ilm and as - The error 
bars on the points indicate the 68% confidence regions. The bias is 
no longer statistically significant if the posterior probability for wq 
peaks well within one sigma from the reference value wq = — 1. 
For the AGN model this is only achieved when s < 500 (this cor- 
responds to a real space separation of more than 40 arcminutes!). 
However, drastically limiting the range of scales increases the sta- 
tistical uncertainty by almost a factor of three. We therefore con- 
clude that this approach is not viable and that one needs to account 
for the effects of baryon physics when computing the constraints 
on cosmological parameters. 



5 REDUCING THE BIAS USING A SIMPLE MODEL 

The results presented in the previous section suggest that one can- 
not ignore the effects of baryon physics on the matter power spec- 
trum in the case of future lensing surveys. A complication is that the 
bias itself depends strongly on the details of the feedback model, 
but that we do not know for sure which of the feedback scenar- 
ios (and parameters) is correct. However, baryon physics also has 
an impact on other observables, which can be used to discriminate 
between models. 

For example, McCarthy et al. (2010) showed that the AGN and 
REF simulations yield haloes with significantly different gas frac- 
tions. Similarly, the amount of gas that cools to form stars is dif- 
ferent, leading to different luminosities. For both observables, the 
AGN simulation provided a good match to observations of groups 
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of galaxies whereas the REF simulation did not. In principle, such 
observations can be used to select hydrodynamic simulations that 
best describe our Universe. In this section, we will explore a dif- 
ferent approach and show that those same observables can be used 
to modify the dark matter power spectrum such that it accounts for 
most of the effects of baryon physics. 



5.1 Halo model 

To predict the matter power spectrum analytically, we take advan- 
tage of the fact that the clustering of haloes of a given mass is 
known in the linear regime and that the average density profiles of 
dark matter haloes are specified by their mass. As shown by Seljak 
(2000), this "halo model" approach can reproduce the power spec- 
trum into the non-linear regime, although some parameters have to 
be calibrated using numerical simulations. 

The power spectrum is computed as the sum of two terms. The 
first one describes the correlation of the density fluctuations within 
the same halo. This Poisson term P^{k) dominates on small scales 
and is given by 



(27r)3 



M{v) 



y[k,M{u)]\ (12) 



where p is the mean matter density and y[k, M{v)] is the Fourier 
transform of the density profile of a halo with virial mass M{v) 
normalised such that: 



y[k,M] = 



JJ"^" 4:nr'^drp{r) 



(13) 



where p(r) is the density profile of the halo and r^ii its virial radius. 
The peak height u of such an overdensity is defined as 



Scjz) 
a{M) 



(14) 



where 5c is the linear theory value of a spherical overdensity which 
collapses at a redshift 2. a{M) is the rms fluctuation in spheres that 
contain mass M at an initial time, extrapolated to z using linear 
theory. We use jSheth & Tormen 1999> : 



(15) 



where u' = af with a — 0.707 and p — 0.3. The normalisation 
constant A is determined by imposing J f{i')di' — 1. Note that the 
function /(u) is related to the halo mass function dn/dM through 



^dM = ■^f{u)du. 
dM Ar^ ' 



(16) 



The second term, P^^{k), describes the clustering of haloes and 
dominates on large scales. It is given by 



P'^'^W = Piin(fc) / duf{u)h{p)y[k,M{y)] 



(17) 



where Piin(fc) is the linear power spectrum, and the halo bias h{v) 
is given by (Mandelbaum et al. 2005): 



h{u) = 1 



u' - 1 



2p 



(18) 



with a = 0.73 and p — 0.15. Finally, one generally assumes that 
the density profile is of the form ( INavarro et al. 19931 : 



p(r) oc 



1 



(19) 



where is the scale radius. Numerical cold dark matter sim- 
ulations have shown that this NFW profile is a fair description 
of the radial matter distribution for haloes with a wide range in 
mass. They also indicate that is not a free parameter, but that 
it is related to the virial mass (albeit with considerable scatter). It 
is customary to account for this correlation by specifying a rela- 
tion, between the concentration c = Vvh/rs and the virial mass. 
We use the mass-concentration relation derived by Duffy et al. 
(2008). The virial mass and radius are related through Alvh- ~ 
(47r/3)r^ir/5c<5vir, where we use the fitting formula of Bryan & 
Norman (1998) to compute Svii{z). 

It is typically assumed that the total matter density can be de- 
scribed by the NFW profile. However, if the stars and the gas do not 
follow the dark matter profile, then the resulting mass profile and 
thus the power spectrum will be different. In the remainder of this 
section we will explore whether it is possible to simply modify the 
density profile to better describe the distribution of the baryons, and 
whether this model can reduce the biases discussed in the previous 
section. Because of the back-reaction of the baryons on the dark 
matter, baryonic effects will also induce changes in the dark matter 
distribution within haloes (e.g., Duffy et al. 2010) and in the dark 
matter power spectrum (van Daalen et al. 201 1). For simplicity, we 
will, however, ignore this complication. 



5.2 An improved halo model 

The halo model is only a good description of the total matter distri- 
bution if the baryons trace the underlying dark matter distribution. 

As feedback processes redistribute the baryons (which make 
up ~ 17% of the total amount of matter), we expect the total matter 
power spectrum to be modified. In this section we explore whether 
it is possible to predict the power spectrum from hydrodynamic 
simulations using 'observations' of the gas fraction and stellar mass 
(or luminosity). 

Somogyi & Smith (2009) showed that in order to construct ac- 
curate models that include CDM and baryons one should treat the 
density field as a two-component fluid. The procedure that we use 
here, which is to compute the linear matter power spectrum using 
the Eisenstein & Hu (1998) transfer function, corresponds to a one- 
component fluid approximation. Somogyi & Smith (2009) showed 
this approximation leads to a biased power spectrum, in particu- 
lar, for high redshifts and for the baryonic component. However, 
following their studies, the bias is less than 1% for z < 3. Since 
this effect is much smaller than the effect we are interested in, we 
can neglect it. Furthermore, we assume that any modification of 
the final power spectrum is caused by phenomena which affect the 
baryonic component in haloes and happen after the halo has col- 
lapsed. This implies that haloes collapse in the same way as in a 
CDM-only Universe. We know this assumption is not completely 
correct, as van Daalen et al. (201 1) have demonstrated that the back 
reaction of the baryons onto the CDM causes the power spectrum of 
the CDM component to differ from that of a dark matter only sim- 
ulation. However, the effect of the back reaction on the total matter 
power spectrum is smaller than the effect caused by the change of 
the distribution of baryons. 

We start by modifying the density profile to better describe 
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the distribution of tlie stars and the gas. Note that the NFW profile 
with the mass-concentration relations of Duffy et al. (2008) is still 
used to describe the dark matter distribution. The stellar mass in a 
halo is much more concentrated than either the gas and dark mat- 
ter components, and we therefore approximate its distribution by 
a point mass. To describe the gas component, we use a single /3- 
model, which provides a fair description to X-ray observations in 
groups and clusters of galaxies (e.g., Cavaliere & Fusco-Femiano 
1976; Reiprich & Bohringer 2002; Osmond & Ponmann 2004). The 
corresponding density profile pgas, is given by 



L V arson/ J 



-3/3/2 



(20) 



where a is the ratio between the characteristic scale of the gas pro- 
file (i.e. the core radius) and rsoo, the radius of a sphere with aver- 
age density 500 times the critical density. The value for the slope 
P is usually considered a free parameter in the fit to X-ray data, 
but we note that for a hydrostatic isothermal sphere its value corre- 
sponds to the ratio of the specific energy in galaxies to the specific 
energy in the hot gas (e.g., King 1972; Jones & Foreman 1984; 
see also Mulchaey 2000 for a review). Finally, note that although 
different models have been used in the literature (e.g., Osmond & 
Ponmann 2004; Amaud et al. 2010), the differences between these 
models are too small to be important here. 

Our model provides a convenient description of the final mass 
distribution because even though the parameters depend on the de- 
tails of the various feedback processes, they can be constrained ob- 
servationally. The most important observable is the gas fraction as 
a function of halo mass: the gas has an extended profile and the 
overall profile depends on the fraction of gas which is left in the 
halo. This has been measured by McCarthy et al. (2010; 2011) for 
the OWLS simulations and effectively determines the value of po in 
Eqn.l20l The results are presented in the left panel of Figure [8j it is 
clear that different feedback models lead to rather different scaling 
relations. 

The values of a and /? have been measured from the simula- 
tions. We note that these values have a large dispersion and depend 
on the mass; however, we take here an average value measured 
over a range of masses 10^^ Mq < A/500 < 10^^ /i^^A/©. 
This corresponds to a = 0.01 and /3 — 0.4 for the AGN and 
DBLIMFV1618 models and a = 0.2 and /3 = 0.9 for the REF 
model. The power spectrum is mildly sensitive to the choice of the 
slope /3 as this slope defines how fast the gas power spectrum de- 
clines for large k. The AGN model, for which a large fraction of 
the gas has been ejected beyond rsoo, has a power spectrum whose 
shape is less sensitive to the choice of f}. 

The right panel shows the relation between A/soo, the mass 
of a sphere with average density 500 times the critical density 
of the universe, and Lks.bcg, the stellar luminosity of the (cen- 
tral) brightest galaxy of the halo in the A'— band as predicted by 
McCarthy et al. (2010) using the metallicity-dependent population 
synthesis model of Bruzual & Chariot (2003). To constrain the pa- 
rameters of our extended halo model, we fit a power law to the data 
in Figure[8]and we find: 

/gas//6(rsoo) = 0.151ogio(Msoo)-1.49 REF 

/gas//6(rsoo) = 0.30 logio(Msoo) -3.40 DBLIMFV1618 

/gas//6(rsoo) = 0.40 logio(Msoo) -4.94 AGN 

logio(ix.,BCG) = 0.621ogio(M5oo) + 4.19 REF 

logio(Lx.,BCG) = 0.811ogio(Msoo) + 1.39 DBLIMFV1618 

logio(Lif.,BCG) = 0.83 logio(Msoo) + 0.23 AGN 



To compute the stellar mass from the stellar luminosity, we measure 
the average mass-to-light ratio of the BCGs directly from the simu- 
lations. We find AL/L^, = 0.65Mq/Kq for AGN, O.32A'/Q/iC0 
forDBLIMV1618, and 0.33Mq/Kq for REF We could have mea- 
sured the stellar mass of the BCG directly from the simulations, but 
since we want to use quantities which can in principle be observed, 
we derive the stellar mass using the same procedure that one needs 
to use for real data. 

To determine the fraction of baryons (stars -1- gas) that remains 
within rsoo, we use the scaling relations shown in Figure |8] We 
then assume that the gas fraction within rvir is the same as within 
rsoo, although we note that may not be correct. The resulting halo 
density profile is the sum of the stellar, gas and dark matter pro- 
files weighted by their respective fractions. We assume that the 
ejected gas is still associated with the the parent halo and make 
ad-hoc assumption that the gas is distributed uniformly between 
rvir < r < 2rvir. We add this unifom component to the function 
y[k, M{h')] so that the total mass of the halo is unchanged. The 
density distribution of the gas at large radii can been derived, using 
X-ray observations, only for bright clusters (Reiprich et al. 2008[ >, 
and it might be overestimated ( S imionescu et al. 201 iK Stacked X- 
ray (e.g., Dai et al. 2010) and SZ observations (e.g., Afshordi et al. 
2007; Komatsu et al. 201 1) could allow one to obtain a more accu- 
rate estimation of the gas distribution even for less massive haloes. 
We will measure the distribution of this gas using the simulations 
and evaluate the importance of its modelling in a future work. 

Under these assumptions our model has more power than the 
simulations for high k. This is due to an overestimate of the stellar 
mass which can be explained as follows. Using the simulations data 
by McCarthy et al. (2010) we can only establish a relation between 
A/500 and Lks.bcg for haloes with A/500 > 5 x 10^^ Mq. 
Observationally the relation between the fraction of stars and halo 
mass is approximately linear down to a certain halo mass; for 
haloes with mass below a few times 10^^ h^^ Mq this relation is 
poorly constrained. However, observations suggests a stellar frac- 
tion for Milky Way size haloes which is about 0.01, suggesting 
that the stellar-mass to halo-mass relation has to change. This ob- 
served relation between luminosity and halo mass agrees well (see 
e.g. Moster et al. 2010) with models constructed using halo occupa- 
tion distributions ( e.g.. Peacock & Smith 2000; Seljak 2000; White 
2001; Berlind & Weinberg 2002) and conditional luminosity func- 
tions (Yang et al. 2003; van den Bosch et al. 2007; Cacciato et al. 
2009), as well as with semi-analytic models (e.g. Croton et al. 2006; 
Bower et al. 2006). Those approaches all show that there is a char- 
acteristic mass A/ ~ 10^^'^ h^^ AIq for which star formation is 
particularly efficient, whereas for small masses the fraction of stars 
in haloes decreases steeply. This should also happen in our simula- 
tions although we lack the resolution required to measure the stellar 
mass-to-light ratio for small haloes. Thus, we just assume that the 
mass of the stars is 0.01 A/500 for A/soo < 5x 10^^ h^^ A/q which 
approximately corresponds to the fraction of stars in a Milky Way 
size halo. Note that while low-mass haloes do not contribute much 
at small k, they will dominate the signal if the stellar fraction is too 
large. 

In Figure |9] we compare the ratios of the power spectra mea- 
sured from the simulations to the ones predicted by the model. De- 
spite its simplicity, the model is able to reproduce the main features 
of the power spectrum for all scenarios for redshifts to which the 
weak lensing signal is most sensitive. Although it is a clear im- 
provement over the pure dark matter model, the model is still too 
simplistic and does not perform equally well at all redshifts and all 
scales. The accuracy of the toy model can clearly be improved and 
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Figure 8. Left panel: relation between A/500 and gas fraction in haloes for REF, DBLIMFV1618 and AGN. Right panel: relation between M500 and the 
stellar luminosity in the K-hand Lxa of the brightest galaxy in the halo. The solid lines represent the best linear fits to the points (see text). The simulation 
data was taken from McCarthy et al. (2010), who showed that the AGN model agrees with the available observations of groups. 




Figure 9. Comparison of the ratio of the power spectra REF/DMONLY, DBLIMFV1618/DMONLY and AGN/DMONLY measured from the simulations 
(solid lines) and predicted by our simple halo model (dotted lines) which has been tuned to fit the gas and stellar mass fractions predicted in each scenario (see 
text). The model matches the simulations fairly well for k < lOh Mpc~^. 
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Figure 10. Contribution of haloes of different masses to tlie one lialo-teim 
{k) (solid lines) redshift z = 0.5. We show for comparison the halo- 
halo term P^^(k) (dotted lines) for the same mass ranges. At scales where 
the baryonic feedback is important, the dominant contribution is the Poisson 
contribution of haloes with M^m) > 10^^ Mq. 



we will aim to do so in future work. In particular, we will examine 
requirements on the types of observations that will be needed to 
constrain the model and the mass and length scales that need to be 
studied. 

Indeed, we remind the reader, that we currently calibrate the 
model using scaling relations covering a large range of masses; in 
particular down to low-mass haloes. Figure [TO] shows the contri- 
butions to the power spectrum from various ranges in halo mass. 
This demonstrates that haloes with masses Msoo > IO^^/i-'Mq, 
i.e. massive galaxies, groups and clusters, contribute most on the 
scales of interest, but that the properties of halos with masses as 
low as 10^^ Mq need to be understood as well. 

Currently, scaling relations can be measured observationally 
down to virial masses of a few 10^^ h-^ Mq, albeit only for bi- 
ased samples. For those high-mass haloes the parameters of our 
model can be constrained directly by observations of groups and 
clusters of galaxies, most notably the gas fraction (e.g., Giodini 
et al. 2009; Sun et al. 2009). Although studies of ensemble aver- 
aged properties can also provide useful constraints (Leauthaud et 
al. 2010; Rykoff et al. 2008), for masses down to Afvir ~ 10^^ Af© 
observational constraints are expected to remain limited. Further- 
more, in our model, we made a rather ad-hoc assumption that the 
gas that is removed from the central regions is in the outskirts of the 
haloes. Although this seems plausible (see for example Rasheed et 
al. 2010) it is not clear how accurately one can observe the distri- 
bution of the gas in the outskirts of the haloes. Because of those 
observational limitations, we expect that an accurate description of 
the baryonic feedback will have to rely in part on simulations. 

Finally, our findings suggests that a model with few parame- 
ters, similar to the one we introduced, can capture the main features 
of the baryonic feedback. Those parameters need to be varied to- 
gether with the cosmological parameters when one performs a cos- 
mological interpretation. One can use observations and simulations 



to derive priors on those parameters and increase the accuracy of 
the cosmological parameters estimation. 



5.3 Improved constraints on cosmological parameters 

To examine whether the biases on the cosmological parameters are 
indeed reduced, we repeat the likelihood analysis described in Sec- 
tion |4l using the improved halo model to compute the lensing sig- 
nal. The results are presented in Figure [TT] For the AGN scenario, 
where the biases in f^m and uiq were ~ 20% and 40% respectively, 
the bias is now visible reduced, although it is still larger than the 
statistical error for a Euclid-like survey. 

Despite its simplicity, the model allows for a significant im- 
provement and its accuracy can be further increased. For instance, 
we did not include any evolution in the scaling relations. We also 
did not account for the intrinsic dispersion in the model parame- 
ters. Furthermore, for the gas that has been ejected beyond rsoo, 
we simply assumed an homogeneous distribution, but more work 
is needed to validate (or improve) this assumption. We will explore 
various improvements to our model in future work, as they are be- 
yond the scope of this paper. 



6 CONSEQUENCES FOR OTHER PARAMETERS 

In section |4] we showed that ignoring the effects of baryons on the 
power spectrum leads to significant biases in the recovered dark 
energy equation of state wq, as well as in cts and f2m. The mea- 
surement of these parameters is, however, not the only goal of fu- 
ture cosmic shear studies. In this section we examine the impact of 
baryon physics on the accuracy of some cosmological parameters 
for which weak lensing has been suggested to be a useful probe. 
These are the scale dependence of the spectral index of the pri- 
mordial power spectrum (section 16. U . the properties of cold dark 
matter (section |6^ and the neutrino masses (section|63}- 



6.1 Running of the spectral index 

The large-scale structure we observe today is believed to arise from 
quantum fluctuations that grew to cosmological scales during a pe- 
riod of inflation. Simple models of inflation predict a primordial 
power spectrum that is nearly scale invariant, P{k) oc fc"^, with 
scalar spectral index ~ 1. In more complicated inflation the- 
ories Ua may also depend on k. A general expansion adopted for 
ns{k) is l |Kosowsky & Turner 1995[ l: 



dink 



In 



(21) 



where fco — 0.015 Mpc^^ is a pivot. A non-zero value for 
Us = dus/d In k results in a running spectral index. Measurements 
of both Us and as provide direct constraints on the inflationary 
models, and are thus of great interest. CMB temperature fluctua- 
tions as well as CMB polarisation measurements are able to pro- 
vide constraints on the spectral index. The most recent results from 
ACT WMAP7 suggest Us (fco) to be slightly smaller than unity and 
Us <0 (Dunkley et a l. 2010| l. 

Weak lensing observations enable us to extend the measure- 
ments to smaller scales, thus tightening the constraints on as. We 
are concerned, however, that such a measurement will also be af- 
fected by baryon physics, because the suppression in the power 
spectrum seen for the AGN model could also be achieved by a 
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Figure 11. Left panel: joint constraints on f2m and erg for a Euclid-like survey, using our improved halo model tuned to fit the gas and stellar mass fractions 
predicted by each scenario to compute the likelihood. The contours indicate the 68% confidence regions for the vaiious feedback models. Middle panel: joint 
constraints on and wq. Right panel: posterior probability distribution of uiq after marginalisation over Cm and erg. When the improved model is used to 
compute the lensing signal, the biases in the recovered cosmological parameters are significantly reduced. 



more negative value for as . This is demonstrated by Figure [T2l 
where we show that it is possible to find dark matter only models 
(i.e. halofit models) with a running spectral index that match the 
simulated power spectra for the lowest redshift bin. 

We also computed the joint likelihood for Ua £ [0.90, 1.20] 
and as £ [—0.10,0.050], marginalising over a flat prior erg G 
[0.65, 0.83] while keeping wq and Q.m fixed. The results are pre- 
sented in the left panel of Figure[T3] As expected, the best-fit value 
for as is biased to a more negative value in the case of the AGN 
model if the effects of baryons are ignored. Our simple extension of 
the halo model is able to reduce the bias also for these parameters, 
as is demonstrated by the right panel in Figure [TJI 



6.2 Warm dark matter 

A number of observations appear to disagree with the outcome of 
ACDM simulations which predict a large amount of structure in 
galaxy-sized haloes. In particular, the simulations predict a larger 
number of satellite galaxies for the Milky Way than is observed (the 
so-called missing satellites problem), although the recent discovery 
of faint satellites appears to alleviate this problem (see Belokurov 
et al. 2010). Another longstanding issue concerns the inner slope 
of the density profile: low surface brightness galaxies show profiles 
which are not as cuspy as the NFW profile (e.g. Kuzio De Naray et 
al. 2010; Oh et al. 2010). 

The existence of warm dark matter (WDM) can, in principle, 
solve the problems mentioned above, because the free- streaming of 
the particles suppresses small-scale fluctuations and leads to haloes 
with a core. We note, however, that Villaescusa Sl Dalai (201 1) re- 
cently showed that WDM haloes cannot form cores that are large 
enough to match the observations. Recent papers suggest that fu- 
ture weak lensing missions such as Euclid, combined with Planck 
observations, can constrain the mass of the WDM particles, through 
the effect they have on the growth of structure (e.g., Markovic et al. 
2010). 

To ensure that structure formation is only altered on small 
scales, the parameters of WDM models are generally chosen such 
that that the power spectrum is modified at subgalactic scales 
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Figure 12. Ratio of the weak lensing power spectra for the lowest redshift 
tomographic bin, for the REF/DMONLY, DBLIMFV1618/DMONLY and 
AGN/DMONLY models (solid lines). The dot-dashed lines show the ratio 
for the dark matter only halofit model that best fits each of the simulations, 
obtained by varying erg, the spectral index ns(fco) and the running as = 
dus/dlnk. In general it is possible to find a combination of parameters, 
such that the dark matter only model with a running of the specfi'al index 
fits the simulated power spectra reasonably well. 



fcfs ~ 0.1 — l/i Mpc~^, which are also affected by baryon physics. 
As an example, we consider the case of fully thermalised particles 
described by Markovic et al. (2010). In this case the WDM streams 
freely shortly before the Universe enters the matter-dominated era. 
The effect on the power spectrum can be parametrised by (Viel et 
al. 2005): 



T{k)={^^y'"=[i+{akf-^r^^-\ 



(22) 



-PcDM(fc) . 

with fi — 1.12 and where PcDM(fc) and PvvDM(fc) are the CDM 
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Figure 13. Left panel: joint constraints on tlie primordial spectral index ns{ka) and running cis = dns/dlnk after marginalising over trg (while keeping 
the other parameters fixed) for a Euclid-like survey. The contours indicate the 68% confidence regions for the feedback models that are indicated. Right panel: 
constraints when we use our extension of the halo model (see Section|5), which was tuned to fit the halo gas fractions and /f-band magnitudes predicted by 
the simulations. The improved halo model partially reduces the level of the bias. 



and WDM power spectra, respectively. The value for the break 
scale a is then given by (Viel et al. 2005): 



a = 0.049 



mwDM 

keV 



ii/fiv 



V 0.25 



r N 0.11 / h \ 1-22 1 

-) Or) '^^^^^ 

(23) 

where f^wDM is the WDM density parameter and ttiwdm is the 
WDM particles mass. Observations of the power spectrum of 
the Lyman-a forest suggest mwDM > 4keV ( |Seljak et al. 2006| 
IViel et al. 2008l|Boyarsky et al. 2009) . If we use this limit, we find 
that for k — Wh Mpc~ , which is where the suppression of the to- 
tal power by baryons is maximal (see Fig.[T}, the power spectrum is 
suppressed by only 5%, which is much smaller than the suppression 
predicted by the AGN feedback model. The suppression decreases 
rapidly, being less than a percent for k — Ih Mpc^^. One could 
argue that this result, which has been derived in the linear approxi- 
mation, is not valid at small scales where the power spectrum is in 
the non-linear regime. However, Boehm et al. (2005), showed that 
the non-linear evolution reduces the difference between CDM and 
WDM models (also see Markovic et al. 201 1). 

We therefore conclude that the effect of WDM on the weak 
lensing signal is likely to be dwarfed by baryon physics. If we wish 
to improve upon current constraints on the WDM mass, we would 
need to be able to compute the effects of baryon physics on the total 
matter power spectrum to less than a percent for k ^ 1 /iMpc^^. 
It is not clear that such accuracy can be achieved. 



6.3 Massive neutrinos 

Particle physics experiments have conclusively shown that neutri- 
nos have mass and suggest that there are three species. Experiments 
are also able to measure the differences in the square of the masses, 
and thus provide a lower bound on the mass of the heaviest neutrino 
species. The effect of massive neutrinos on the power spectrum is 
similar to that of warm dark matter. The main difference is that 
neutrinos are much lighter (i.e, they are considered hot dark mat- 



ter). As a consequence, they become non-relativistic in the matter 
dominated era. 

While experiments are designed to measure the mass of 
the electron neutrino, studies of large-scale structure constrain 
the total density of neutrinos in the universe. The most recent 
estimate from WMAP7 is < 1.3 eV, or "^rrii, < 

0.65 eV if constraints from baryon acoustic oscillations are in- 
cluded ( Komats u et al. 20TT] |. This yields an upper bound for on 
the total neutrino density of fi^ < 0.02 for WMAP7, where we 
used the relation between the neutrino mass and their cosmological 
density: 



93.8/i2eV 



(24) 



As before, massive neutrinos do not modify the power spec- 
trum on scales larger than their free-streaming length, which de- 
creases with time. Therefore there is a fcnr below which the power 
spectrum is not affected ( [Lesgourgues & Pastor 2006) : 



: 0.018 .-(^) 



1/2 



(25) 



To compute the effect of neutrinos on the power spectrum, we use 
the transfer functions given by the Boltzmann code CAMB (Lewis 
et al. 2000). Figure [14] shows the ratio of power spectra at z — 0.5 
with neutrino fractions = = 0.01 and fi, = 0.05. Note 

that the best-fit WMAP7 cosmology JKomatsu et al. 201 iT) yields 
fcnr < 0.01 h Mpc-^ and < 0.1. 

At intermediate scales our results indicate that neutrinos cause 
a scale-dependent suppression of the power spectrum. For large k, 
our calculations, which are based on linear perturbation theory, in- 
dicate that the ratio Pv{k) / PcDM{k) is (almost) constant with a 
value ~ (1-8/^) dEisenstein & Hu 1998l l. A concern might be that 
linear theory does not provide an adequate description. Although 
the highly non-linear regime has not been studied, the numerical 
simulations presented in Brandbyge et al. (2008) suggest that our 
results are sufficiently accurate: Pv{k) / PcDM{k) reaches a mini- 
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Figure 14. Ratio of the AGN/DMONLY power spectra (blue line), and dark 
matter power spectra with = Q^/Qni = 0.01 and 0.05, which cor- 
respond to neutrino masses of '^rni^ ^ 6.0 and "^rrii, ~ 1.2 eV, re- 
spectively. The effect of massive neutrinos on the power spectrum is quite 
different from that of baryon physics, even if neutrinos are light. 

mum value of ~ ( 1 — 10/^ ) and slightly increases again for smaller 
k. The results presented in Figure [14] indicate that for ~ 0.05 
(or a mass of J]] rriu ~ 6.0 eV) the suppression of power by mas- 
sive neutrinos dominates over the baryon physics. In the case of 
light neutrinos (^mi, ~ 1.2 eV or fi, = 0.01), baryon physics 
dominates for fc > 1 /iMpc~^, but the neutrinos cause a signifi- 
cant suppression of power on intermediate scales. Thus, it seems 
likely that one would be able to get useful constraints using only 
large-scale measurements. 



7 CONCLUSIONS 

Weak gravitational lensing by large-scale structure is a powerful 
technique to study the properties of dark matter and dark energy. 
However, the interpretation of the signal requires comparison to a 
predicted matter power spectrum. To date, such cosmic shear mea- 
surements have used predictions based on numerical simulations 
of cold dark matter. As the precision of future projects improves, 
the effects of baryon physics need to be quantified. Recently, van 
Daalen et al. (2011) have shown that the impact of baryonic feed- 
back is significant and demonstrated that the precision with which 
one is able to predict the matter power spectrum from simula- 
tions depends on how well one understands and models mecha- 
nisms such as radiative cooling, star formation, metal enrichment, 
supemovae and AGN feedback. In particular, they have shown that 
baryon physics modifies the matter power spectrum to a level which 
is higher than the accuracy required for future weak lensing exper- 
iments, which require predictions for the amplitude of the power 
spectrum with an accuracy of a few percent for a large range of 
scales. 

In this paper we made a first attempt to quantify the effect 
of baryon physics on the cosmic shear signal, using a subset of 
the cosmological hydrodynamic simulations that were carried out 
as part of the OWLS project (Schaye et al. 2010) and for which 
the matter power spectra have been measured by van Daalen et al 
(2011). 



We find that if weak lensing measurements are interpreted us- 
ing a power spectrum that ignores the effects of baryon physics, 
the resulting cosmological parameter estimates are biased. The bias 
depends on the feedback model, but is as large as ~ 40% for the 
dark energy equation of state wo for the AGN scenario, which we 
consider our most realistic model as it reproduced the observed X- 
ray and optical properties of groups of galaxies (McCarthy et al. 
2010). The biases dominate over the statistical errors for Euclid, a 
proposed space mission that aims to measure wa with a precision 
of a few percent. 

The reason for the large bias is the fact that efficient feedback 
modifies the power spectrum out to several megaparsecs, contrary 
to the commonly held belief that only very small scales are affected, 
as is the case in previous models that suffered from the "overcool- 
ing" problem. To reduce the bias to acceptable levels, one would 
need to discard scales smaller than ~ 40 arcminutes from the weak 
lensing analysis in the case of Euclid, at the expense of a factor 
three increase in the statistical error. 

Such a large reduction in sensitivity to the dark energy prop- 
erties is not desirable. Instead, we explored an approach to reduce 
the bias, based on a simple extension of the halo model. We treat 
the dark matter, gas and stars as separate components, assigning 
suitable profiles to each. The baryonic mass assigned to each com- 
ponent is determined using scaling relations between the halo mass 
and the gas fraction and stellar mass. For each scenario the model 
parameters were constrained by fitting the simulated gas fraction 
and stellar masses of groups and clusters at redshift 2 = 0. We also 
made the ad-hoc assumption that the ejected gas is distributed uni- 
formly between the virial radius and 2rvir. Further work is needed 
to improve upon this assumption. 

Although several improvements should be made (an issue that 
we will explore in future work), the simple model presented here 
is able to describe the OWLS power spectra fairly well for k < 
lO/iMpc"^. Particularly encouraging is the reduction in the bias 
in Wo and the other cosmological parameters. 

The usefulness of our model hinges partly on our ability to 
calibrate the model though observations of the distribution of gas 
and stars in halos with masses as low as 10^^ Mq, which may 
be too low for X-ray and SZ observations. Moreover, in order to 
achieve the precision required for future weak lensing experiments, 
one would need to know the distribution of the gas in the peripheral 
regions of haloes. This suggests that improvements require a 'hy- 
brid' approach, relying on a combination of deeper observations 
and better simulations. For instance, dedicated X-ray and SZ stud- 
ies can help calibrate the model parameters, but they also provide 
constraints on the accuracy of the hydrodynamic simulations. 

Finally, we also examined the effects of baryon physics on 
other applications of weak lensing. For Euclid, constraints on the 
parameters of the primordial power spectrum predicted by theories 
of inflation are also significantly biased. The measurement of the 
neutrino mass is also affected (in particular if the mass is small), 
but the bias should be small if only scales larger than ~ 10 Mpc are 
considered. Improved weak lensing constraints on the properties of 
warm dark matter are unlikely, because current constraints lead to 
changes in the power spectrum that are much smaller than those 
introduced by baryon physics. 
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